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We investigate the muonic hydrogen 2Pj^J^ to 2S^J^^ transition through a precise, non- 
perturbative numerical solution of the Dirac equation including the finite-size Coulomb force and 
finite size vacuum polarization. The results are compared with earlier perturbative calculations of 
(primarily) Borie, Martynenko, and Pachucki 0-01; and experimental results recently presented by 
Pohl et al. 01 , in which this very comparison is interpreted as requiring a modification of the proton 
charge radius from that obtained in electron scattering and electronic hydrogen analyses. We find 
no significant discrepancy between the perturbative and non-perturbative calculations, and present 
our results as confirmation of the perturbative methods. 

PACS numbers: 36.10.Ec,31.30.jr,03.65.Pm,32.10.Fn 
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I. INTRODUCTION 



The precision measurement of the Lamb shift tran- 
sition energy between the 2P^^^ and 2S^j^^ states of 

muonic hydrogen by Pohl et al. 0, see Fig. [TJ has cre- 
ated considerable interest because of a 0.31 meV discrep- 
ancy with the value predicted by theoretical calculations 
(speficically those discussed in Ref. 0; Borie H-Q and 
Martynenko d Q along with many others j8l4lll|). This 
Lamb shift splitting of 0(206) meV is dominated by the 
lowest order QED vacuum polarization, and obtains a 
significant contribution from the finite size of the proton. 
On selecting and combining the perturbative predictions 
for the corresponding contributions to the measured tran- 
sition, Pohl et al. produce a cubic equation relating their 
experimentally measured energy shift to the theoretical 
prediction, and arrive at 

206.2949(32) meV = 206.0573(45) 

-5.2262(r2)i/2 
-F0.0347(r2)3/2 meV, (1) 

the only physically-meaningful solution of which im- 
plies a proton rms charge radius of = ^ (r^) ~ 

0.84184(67) fm which differs from the consensus CO- 
DATA [ill value of rp = 0.8768(69) fm by 4.9 standard 
deviations. 

Such a large modification of a basic electromagnetic 
property of the proton suggests that either there may 
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be an as yet unrecognised problem in several other ex- 
perimental efforts (such as the electronic hydrogen spec- 
troscopy and scattering experiments which primarily lead 
to the CODATA value [l3) or in the QED calcula- 
tions [isl - fTsj , or alternatively that some new physics (be- 
yond the Standard Model) contributes to this transition 
energy [l6l - [T8j . With respect to the QED calculations, 
we note that the predominant theoretical approach in- 
volves perturbation theory applied to the solutions of the 
non- relativistic Schrodinger equation [H, 0, Q • Since the 
effects of finite size and the vacuum polarization potential 
are quite large at short distance, it seems important to 
verify by explicit calculation that the perturbative treat- 
ment is indeed adequate at the level quoted. 

We therefore calculate the transition energy relevant 
to the aforementioned experiment using the relativistic 
Dirac equation to describe the muon wave function non- 
perturbatively. We take care to control the numerical 
errors in the calculations and to quantify any differences 
from the perturbative non-relativistic approach. In our 
work we extend considerably the earlier work by Borie P, 
i- 

In the sections following, we discuss the nature of the 
transition and contributing physical effects, as well as the 
method by which we calculate the energies corresponding 
to the various eigenstates. We summarize any discrepan- 
cies with respect to previous work. 



II. MUONIC HYDROGEN SPECTRUM 

After muon capture by hydrogen, about 1% of the 
muons reach the metastable 2S'-state. This state is acti- 
vated by laser excitation from 25*^2^ ^o 2P^=^ and the 
signature that the laser energy is well tuned is the appear- 
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25'i/2 eigenstate, this leads to hyperfine eigenstatcs with 
total angular momenta F = and F = 1. The muon 
involved in the measured transition decays to the 25^^^^^ 
state, and thus we address an accurate determination of 
the energy of this state and the corresponding splitting 
Similarly, the hyperfine states are labeled by 
and F = 1 for the 2Pi/2 case, and F = 1 and 



F 
F 



2 for the 2P3/2 case. The hyperfine splitting ener- 
gies are denoted here by AE^pg'' and AE^^p'g for the 
corresponding states. The muon in the measured transi- 
tion decays from the 2P^2 ^ eigenstate, so we require an 
accurate determination of the energy of this state. We 
note that the muon-proton tensor force does lead to some 
mixing of the two F=l states but we will not make a new 
calculation of this effect. 



III. NUMERICAL METHOD 



FIG. 1: (Color online) Muonic hydrogen L-shell spectrum, 
showing the finite proton-size correction. Lamb shift splitting, 
fine structure, and hyperfine structure. The 2S^J^^ to 2P^J^ 
transition measured by Pohl et al. [21 is shown as the green 
dotted hne (marked ^). 



ance of a prompt El transition X-Ray of a muon from 
the 2P^^^ state of muonic hydrogen to the IS^^"^^ state. 
Given the resonance method used to establish the energy 
of the Lamb shift, there appears to be no doubt about 
the remarkable experimental result of Pohl et al. . 

The L-shell level scheme explored in this experiment 
is depicted in Fig. [T] for the point-Coulomb potential of 
magnitude 
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The 25i/2 and 2Pi/2 states are degenerate even in the 
Coulomb-Dirac theory. Since the muon orbit is 200 times 
smaller than the corresponding electron orbit, the S- 
states probe the point charge smearing effect of vacuum 
polarization, which is the dominant component in the 
Lamb shift, AP^fmb^' ^^S- E similar nature is a 
smaller but significant effect for hydrogen {Z = 1) arising 
from the finite proton charge radius, denoted ^Epf^^^^ in 
Fig. [TJ While the vacuum polarization effect increases 
the binding of the S1/2 states pulling these 'down' in 
Fig. [U the finite charge distribution acts in the opposite 
direction. 

The spin-orbit fine structure splitting of the 2Pi/2 and 
2P3/2 eigenstates is predicted by the Dirac equation. It is 
denoted by AEpg in Fig. [TJ As the measured transition 
does not involve the 2Pi/2 state, we shall only calculate 
the relevant energy levels for the purpose of comparison 
with perturbative results. 

The spin-spin coupling of the muon to the proton adds 
hyperfine splitting to the spectrum. For the case of the 



To calculate the theoretical energy difference cor- 
responding to the measured transition, previous au- 
thors have primarily used perturbation theory with non- 
rclativistic wave functions, including in the effective in- 
teraction terms describing the various relativistic correc- 
tions. In order to calculate the perturbative effect on the 
energy produced by an operator, SV, additional to the 
Coulomb potential. 
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we require a wave function to integrate over. If we fol- 
low the methods of previous authors, we can use exact 
Schrodinger wave functions for states with quantum num- 
bers n, m and the lowest order correction to the energy 
is 



AP{^f™ 



&Z.ir)SV r/Z.,Ar)d'r. (4) 



An alternate approach, which we choose here, is to 
use the Dirac equation with the appropriate potential in 
order to calculate the perturbed wave functions. This 
approach is known to be a specific limit approximation 
to the two-particle Bethe-Saltpeter equation Q . 

We consider the muon wave function in state a to be 
a spinor ipa(f') 



(5) 



normalised to unity, such that the probability is 

[gc^ir)^ + Urf] dr = 1, (6) 



noting that Xk are eigcnfunctions of the total angular mo- 
mentum operator (consisting of a combination of spheri- 
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cal harmonics and Pauli spinors) satisfying 



mt m' (]f — X ,X 



(7) 



The separation of center of mass motion is not exact 
for a relativistic two body system. To lowest order in the 
ratio of muon to proton mass we use the reduced mass fi 
in place of the muon mass in the Dirac Equation 



Mr, 



(8) 



Along with further recoil corrections treated in perturba- 
tion theory, this should provide a very accurate descrip- 
tion of the system. 

The binding of the muon in this system is extremely 
weak (on the scale of rest energy) and as such the eigen- 
value, Eq, for each state calculated using the Dirac equa- 
tion is approximately equal to the reduced mass /x. In 
order to precisely calculate the variance from this value, 
we rewrite the Dirac equation to incorporate the shift 
of eigenvalue down by the reduced mass, such that the 
eigenvalue we arc now solving for is Aq = — fi, corre- 
sponding to the binding energy. 

The Dirac equation we wish to solve is therefore given 

by 



-K+v — 

r 
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where the value of specific to each eigenstate, 

namely 
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H= -1, 
H = +1, 



25'i/2 : K= -1, 



2F. 



3/2 



K. = -2. 



In order to integrate Eq. ([9]), we supply an initial guess 
for the eigenvalue Aq, and appropriate boundary condi- 
tions for the upper and lower components of the wave 
function at small and large radii. We then integrate 
from each limit towards a central matching point. The 
normalised discontinuity in the wave function integrated 
from each limit is used as a measure of the inaccuracy of 
the eigenvalue and a refined estimate is calculated. This 
process is iterated until Aq, changes by less than the re- 
quired tolerance. 



IV. QUALITY CONTROL 

To convince ourselves that our method is self- 
consistcntly accurate, we check the accuracy of our pro- 
cedure using several methods. The unperturbed point- 
Coulomb Dirac eigenvalues are known analytically to be 



("q - |Kq| + - Z2a2)2 



(10) 



where denotes the principle quantum number for state 
a. We first ensured that were are able to reproduce these 
values within a reasonable compute-time. For the 25*1/2 
wave function, we reproduce the analytic result to within 
10 neV, the 2Pi/2 eigenstate to within 40 neV and the 
l'S'i/2 and 2P3/2 eigenstates to within 10 peV, using quad- 
precision Fortran, a sufhcicntly large grid size, and sufh- 
ciently small grid spacing. 

This test does not assure that solutions for a realis- 
tic case — including the finite size of the proton as well 
as finite size vacuum polarization — converges with same 
accuracy. For this reason we employ the virial theorem 
test for our solutions (refer to Ref. [Ol for further details) 
by calculating the reduced eigenvalue as 

Aq - {ipo^lfiiP - 1) + V{f) + r ■ W(r)I^Q). (11) 

The virial theorem provides a far more stringent test 
of the accuracy of the muon wave function near the ori- 
gin, where \VV\ is greatest. We find that the eigenval- 
ues calculated using Eq. pU)) and Eq. (fTTj) for the 25*1/2 
wave function differ by 180 neV for a point-Coulomb 
potential, and 450 neV the finite-Coulomb plus finite- 
vacuum polarization potentials discussed in Sec. I VII We 
therefore conservatively take our errors to be less than 
±500 neV. Propagating this error, we find that in princi- 
ple we could determine the required proton rms charge- 
radius to within approximately 0.05 am (5 x 10~^ fm). 
This should be sufficient to provide a reliable, indepen- 
dent test of the accuracy of the perturbative approach, 
however we note that the determination of the proton 
rms charge radius cannot be performed to this precision 
as the error in that analysis is dominated by experimental 
error in the determination of the transition energy. 



V. 25i/2-2Pi/2 LAMB SHIFT 

The Lamb shift is the splitting of the otherwise degen- 
erate 25i/2 and 2Pi/2 eigenstates attributed to the vac- 
uum polarization potential Vyp i which for a point source 
is given in Q as 



Wp(r) 
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(12) 

where rrie is the electron mass. We can calculate the 
effect that this has on the eigenvalues by assuming that 
this potential is a small perturbation of the Coulomb 
potential, and thus using Eq. Q we find 



A Tpnlm 
^^Lamb 



VyAr)mZ-oAr)Vd'T, 



(13) 



to which we must also add higher order perturbation the- 
ory, relativistic, recoil, and radiative corrections and gen- 
erally higher-order (in a) corrections. 

Alternatively — and more accurately — we can calcu- 
late the shift in eigenvalues using converged Dirac wave 
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functions in response to the combined effect of the 
Coulomb and vacuum polarization potentials. In this 
case we simply take the difference between the converged 
eigenvalues for the 251/2 a-nd 2Pi/2 eigenstates calculated 
in the presence of point-Coulomb and point-vacuum po- 
larization potentials 



as-2P 

Lamb 



= A 



2Pl/2 
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2Si/2 



= 205.1706(5) mcV, (14) 



Care must be taken when comparing this calculation 
to that of Eq. (fT3|) since our calculation includes rela- 
tivistic corrections, which are treated as corrections to 
Eq. (|13p in Ref. [l| . A summary of this comparison and 
the calculated values for the Lamb shift are given in 
Table HI where we note that the perturbative and non- 
perturbative calculations are found to be in good agree- 
ment. For this table and those that follow, we refer to 
various iterations of our calculations in which we com- 
pute the wave function in the presence of point-Coulomb 
(C); finite (size nucleus)-Coulomb (FC); point vacuum 
polarization (VP); and finite (size nucleus) vacuum po- 
larization (FVP) potentials. The dependence on (ry)" is 
extracted in each case by fitting the energy shifts calcu- 
lated at various values of Vp = {rp)^^^ to a cubic of the 
form given in Eq. for the case of a Coulomb-only 
potential, and with the addition of a term proportional 
to 1 when including the vacuum polarization (to account 
for the 2S-2P splitting). The values listed in the column 
'Pohl et al.' of Table |I] have their origins in references 
1-5, 11, 12, 14-17, and 19-25 of Ref and are taken to 
be reliable values. 



VI. PROTON FINITE-SIZE CORRECTIONS 

The perturbative leading-order contributions associ- 
ated with the finite size of the proton arise from consid- 
eration of the proton-form factor. These are introduced 
by Boric and considered in Friar [l3|, Section VI be- 
low Eq. (64b). These terms appear in Pohl et al. Q as 
quoted from Ref. [l|, and are given by 
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(15) 

To calculate this effect in our fully relativistic, non- 
perturbative calculation, we consider the replacement of 
the point-Coulomb potential with the finite-size Coulomb 
potential in Eq. ([9]) 



, , Za 

Vc{r)^ > 

r 



-Za 



\r — r'\ 
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where p{r) is the proton charge-distribution (or more 
accurately, the Fourier transform of the Sachs electric 
form- factor). In response to concerns that the shape of 



the form factor may significantly influence the theoreti- 
cal calculations [13, we have studied the dependence of 
the finite-size correction on the form of this term (always 
normalised to unity) and this will be summarized in an 
upcoming publication (Ref. [2^ ) , though the dependence 
on the choice of charge-distribution — whether it be ex- 
ponential, Yukawa, or Gaussian in form — appears to be 
extremely weak. 

The exponential form for the charge-distribution, nor- 



malised to unity such that / p{r) d r = 1 is given by 



p{r) = Ae-. 



7? = ^12/ {vD . (17) 



We calculate the Lamb shift by taking the difference 
between the appropriate eigenvalues calculated using the 
Dirac equation with the potential given by Eq. (fT6|) with 
the charge-distribution given by Eq. (|17p for various val- 
ues of (rp). We then interpolate the energy shifts and fit 
the data to a cubic of the form 



(18) 



which provides the relevant parameterization. The (r^)" 
dependence in the presence of an exponential finite-sized 
Coulomb potential and point vacuum polarization poten- 
tial 



V{r) = -Za J ■ 



p(r') 



^dV + Wp(r), (19) 



is given by 

A-Epinito = 205.1706- 5.2169(r2) + 0.0353(r2)3/2 meV. 

(20) 

A further important effect of the finite size of the pro- 
ton arises through the convolution of the vacuum po- 
larization potential (Eq. (|12p ) with the proton charge- 
distribution. This leads to the replacement of the point 
vacuum polarization potential by 



Vypir) 



2Za' 



p(r') 



^Zo(|r-r|)dT', (21) 



where we use the expression given in Ref. [ll| 

^n(k1) = 



1 \ (£2 _ 



and where A denotes the electron Compton wavelength 
(divided by 27r), A = 386.15926459 fm. When discussing 
this potential, it should be assumed that we are using 
a normalised exponential charge-distribution. We once 
again calculate the eigenvalues using various values of 
(rp) in the charge-distribution and fit the resulting ener- 
gies to a cubic (as per Eq. (O), except that in this case 
the vacuum polarization induces the Lamb shift, and we 
must include a term proportional to 1. Thus we find 



Ai^Finitc = 205.182 - 5.2519(r,i) + 0.0546(r 



2\3/2 



meV, 
(23) 
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TABLE I: Contributions to the 2S-2P Lamb shift with comparison to values presented in Pohl et al. 01 which themselves 
are selected values from various theoretical sources — references 1-5, 11, 12, 14-17, and 19-25 of Pohl et al.. Values are all 
in meV. Errors in the Dirac calculations are taken to be ±500 neV as per Section IIVI We refer to various iterations of our 
calculations in which we compute the wave function in the presence of point-Coulomb (C); finite-Coulomb (FC); point vacuum 
polarization (VP); and finite vacuum polarization (FVP) potentials. The dependence on (r^)" is extracted in each case by 
fitting the energy shifts calculated at various values of Vp to a cubic of the form given in Eq. (|18p for the case of a Coulomb-only 
potential, and with the addition of a term proportional to 1 when including the vacuum polarization. The listed corrections 
are already included in our Dirac calculations, namely lines 3 and 5 of Table 1 in Ref. y\ and the nuclear size contributions 
of Table 2 from that reference. All further corrections to both the perturbative calculation and our calculation are contained 
in 'Remaining Corrections' which in this case encompasses all remaining contributions of Table 1 and radiative correction of 
Table 2 of Ref. 0]. 



Contribution 


Pohl et al. 




Present Work 


Dirac (y = Vc -I- Vyp) 
Dirac {V = Vpc) 
Dirac {V = Vfc + Wp) 




205.1706 


205.1706 
-5.2000 (r|> -1- 
- 5.2169 (r|> + 


0.0350 (r-2)3/2 
0.0353 (r|>3/2 


Dirac (V = Vfc + Vpvp) 
Relativistic one loop VP 

Polarization insertion in 
two Coulomb lines 

Finite size effects 


205.0282 
0.1509 

-5.1987 (r|) + 0.0347 (r2)3/2 


205.1822 


- 5.2519 (r|> 4 


- 0.0546 (r|>3/2 


Subtotal: 


205.1791 - 5.1987 (rj) + 0.0347 (r|>3/2 


205.1822 


- 5.2519 (r|) 4 


' 0.0546 (r|>3/2 


Remaining Corrections 


0.8782 - 


- 0.0275 {r2) 






Total: 


206.0573 - 5.2262 (rj) + 0.0347 {rjf/'^ 


206.0604 


- 5.2794 (r|) 4 


' 0.0546 (r|>3/2 



which is the expression vifhich is compared to the per- 
turbative calculation in Table HI We note that the finite- 
vacuum polarization induces a small but non-trivial shift, 
and that the results are otherwise essentially the same as 
those of Pohl et al. 0. 

VII. 2P FINE STRUCTURE 

The 0{Za)'^ perturbative 2P fine structure splitting is 
calculated in Ref. Q to be 

AErs = ^^(l+'-^), (24) 
32ni-^ \ nip J 

along with higher-order corrections. Taking this splitting 
as the difference between the converged eigenvalues of the 
2Pi/2 and 2P3/2 eigenstates gives 

AEjf'g = -^aPj/a - , (25) 

which we can also calculate in the presence of the various 
potentials. For the case of an exponential finite-Coulomb 
potential with finite vacuum polarization, the 2P fine 
structure splitting is 

AEfs = 8.4206(5) meV. (26) 

A comparison of this value with perturbative calcula- 
tions of Borie [J is presented in Tabic HIl The effect of the 



finite-size Coulomb potential (as compared to the point 
case) is negligible at the level of errors of our calculation. 
Similarly, the effect of the finite-size vacuum polariza- 
tion is also negligible at our level of errors. The vac- 
uum polarization itself increases the fine structure split- 
ting by 5 fieY. We note that the perturbative and non- 
perturbative calculations are in perfect agreement to the 
level of errors presented here. 



VIII. HYPERFINE STRUCTURE 

The hyperfine structure is a measure of the £ ■ a cou- 
pling. Following the lead of Ref. [2l|, the appropriate 
Hamiltonian is given by 

n = 2/3,7?i^i^/l\l.J+i|^/37ft|^(0)pI.S, (27) 
j(j-|-l)\H/ 3 

comprising a dipole term and a contact term, for which 
the following definitions apply for the muon Bohr mag- 
neton /3^; proton Bohr magneton /3p; and proton gyro- 
magnetic ratio 7 

/3f, = V^/2to^, Pp ^ ^/^/2Mp, 7 = 2(1 + K)f3p . 

(28) 

Here k ~ 1.792847351 is the proton anomalous magnetic 
moment. V-'(O) represents the muon wave function at the 
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TABLE II: Contributions to tiie 2P fine-structure splitting 
with comparison to values found in Borie Subscripts are 
defined in Table [T] Values are all in meV. Errors in the Dirac 
calculations are taken to be ±500 neV as per Section HVl The 
listed correction (Uehling/vacuum polarization) is already in- 
cluded in our Dirac calculations. All further corrections to 
both the perturbative calculation and our calculation are con- 
tained in 'Remaining Corrections' which are detailed in Ta- 
ble II of Ref. [ll]. Finite-size effects in either the Coulomb 
or vacuum polarization potentials provide no shift above the 
level of errors here, as expected for P-states. The perturbative 
calculation prediction is reproduced within errors. 



Contribution 


Borie Present Work 


Dirac (V = Vc) 


8.4156 8.4156 


Dirac {V = Vpc) 


8.4156 


Dirac (V = Vc -1- Vyp) 


8.4206 


Dirac (V = Vpc + VVp) 


8.4206 


Dirac {V = Vpc + Vpvp) 


8.4206 


Uehling (VP) 


0.0050 


Subtotal 


8.4206 8.4206 


Remaining Corrections 


-0.06852 


Total: 


8.3521 8.3521 



centre of the proton. We now investigate the two terms 
of Eq. (|27)) separately. 



A. 2Si/2 Hyperfine Structure 

There exist several methods by which the 25* hyper- 
fine structure can be calculated. The perturbative 25' 
hyperfine structure calculated in Ref. Q is given by 



AE 



2Si/2 
HFS 



-{Za) 

3 nium^ 



(1 



(29) 



For £ ~ the contact term in the Hamiltonian (Eq. ([27]) ) 
is non-zero, while the dipole term vanishes; 



E^fps = ^(3,lh\m\'{FrnF\I ■ S|Fmj.) , 



(30) 



where iFmp) is the eigenfunction belonging to F = 1 + J, 
such that 



{FmplI ■ SlFmp) 



F(F+ 1) 



(31) 



Thus, the splitting between the 25 F 
hyperfine levels is given by 



and F = 1 



. p2S{F=l-F=0) 



i^/3^7ft|^(0)P, (32) 



We note an important, relevant typographical correc- 
tion; In Ref. [13 Eq. (18.2-17b), the sign should be pos- 
itive and the second 9 in the denominator should not 
appear. 



The value of the 25* hyperfine splitting, as calculated 
using Eq. (|32p with the wave function calculated with 
the Dirac equation in the presence of the combined ex- 
ponential finite-Coulomb and finite vacuum polarization 
potentials is 



AE 



25 
HFS 



22.7690(5) meV. 



(33) 



We note that the effect of including the exponential 
finite-size Coulomb potential as compared to the point 
case reduces the splitting by 0.1269(5) meV; introduc- 
ing the point vacuum polarization potential increases the 
splitting by 0.0747(5) meV for the point-Coulomb, and 
0.0742(5) meV for the finite-Coulomb cases. Using the 
combined finite vacuum polarization and finite-Coulomb 
potentials reduces the splitting by 0.0012(5) meV to give 
the value above. 

Alternatively, one can follow Ref. [H in which case we 
can calculate this splitting to be 



Kgfi 



- [A(A + l)-/(/ + l)-j(j + l)] 



2M„ 



r '^g[r)f{r)dr. 



In that case, we calculate 

^Elifps 22.7640(5) meV, 



(34) 



(35) 



provided we correct for the reduced mass in the formula 
such that the magnetic moment of the muon is not de- 
fined in terms of reduced mass but rather defined in terms 
of the free space mass. 

A comparison to perturbative calculations is given in 
Table Unl where we note the finding of a finite-size depen- 
dent contribution in this splitting, which is neglected in 
the summary of Pohl et al. (though finite-size effects sans 
a parameterization are calculated in studies by Borie [l| 
and Pachucki ) and which differs from results obtained 
using the standard Zemach treatment Q. Wc note that 
for (r2)i/2 = 0.8768 fm, the 251/2 hyperfine splitting 
is calculated here to be 22.8496 meV (22.8547 meV for 
(r2)i/2 ^ 0.84184 fm) which indicates a 0.0087(5) meV 
(0.0100(5) meV) correction to the perturbative calcula- 
tion (once the factor of 1/4 is taken into account), or 
2.8-3.2% of the 0.31 meV quoted discrepancy. 



B. 2Pi/2 Hyperfine Structure 

The 2Pi/2 hyperfine structure is of no consequence for 
the transition which we are investigating here. Nonethe- 
less, we calculate the energy of the 2^^^^" and 2Py'^^ 
levels as a confirmation of our method, and to compare 
to the perturbative results. Following Ref. to C(a^) 
the 2Pi/2 hyperfine structure splitting is given by 



^^HFS 



Ef 



, m^(l-|-2K) 



6 12TOp(l + k) 



(36) 
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TABLE III: Contributions to the 25i/2 hyperfine splitting calculated via Eq. I|32p with comparison to values found in Marty- 
nenko [j. Subscripts are defined in Table HI Values are all in meV. Errors in the Dirac calculations are taken to be ±500 neV 
as per Section HVl The listed corrections are already included in our Dirac calculations and are listed by their descriptions in 
Ref. AH further corrections to both the perturbative calculation and our calculation are contained in 'Remaining Correc- 
tions' which encompasses the muon AMM, amongst other corrections listed in Ref. 4]. We note that the 'Proton structure 
corrections of C(a^)' pertains to the Zemach contribution (which we shall explore in an upcoming publication) and does not 
include considerations of finite-size in the wavefunction, and that the polynomial dependence on (rp)" of this splitting is not 
discussed in the literature. 



Contribution 


Martyncnko 




Present Work 


Dirac (V = Vc) 






22.8229 




Dirac (y = Vc -1- Vyp) 






22.8976 




Dirac {V = Vpc) 




22.7774 


- 0.1746 (r2) -f 


0.0709 


Dirac {V = Vpc -I- Vyp) 




22.8510 


- 0.1701 (r2) + 


0.0667 (r|)3/2 


Dirac {V = Vpc + V'fvp) 




22.8521 


- 0.1795 (rj) + 


0.0739 (r|)^/2 


Fermi Energy Ep 


22.8054 








Relativistic correction ^(Za)^Ep 


0.0026 








VP corrections of orders a^, 










in the second order of perturbation series 


0.0746 








Proton structure corrections of order 


-0.1518 








Proton structure corrections of order 


-0.0017 








Subtotal: 


22.7291 


22.8521 


- 0.1795 (r2) + 


0.0739 (r|)='/2 


Remaining Corrections 






0.0857 




Total: 


22.8148 


22.9378 


- 0.1795 (r2) + 


0.0739 (r|>3/2 
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for which the Fermi energy is 



{ZaY 



(37) 



and where is the muon anomalous magnetic moment. 
We note another important, relevant typographical cor- 
rection; in Ref. [1] the factors of 2 in the denominators 
of the third terms of Eqs. (27-28) should read 12. The 
calculations are performed correctly however. For f 7^ 0, 
the dipole term in the Hamiltonian (Eq. ([27])) is non- 
zero, while the contact term vanishes. The energy for 
the dipole term is thus given by 



' \){FmF\l-3\FmF), (38) 



j{i + 1) \r 



where the non-zero terms in the dot-product are given 
by 

{FmF\l-i\FmF) = i[F(F+l)-/(/+l)-j(j + l)]. (39) 

For Schrodinger wave functions, the vacuum expectation 
value of is analytic, in that 



1 



(40) 



Inserting the appropriate values of F, n, £, I, j for each of 
the F = and F = 1 states, one obtains the energy of 
the 2Pj/2 hypcrfinc structure to be 



(41) 



TABLE IV: Contributions to the 2P1/2 hyperfine splitting 
with comparison to values found in Martynenko 5]. Sub- 
scripts are defined in Tabled Values are all in meV. Errors 
in the Dirac calculations are taken to be ±500 neV as per 
Section lTVl The listed corrections are already included in our 
Dirac calculations. The 'Leading contribution' is extracted 
from the ©(Za)* line in Table II of Ref. [5[ and includes 
only the Ef/3 term of Eq. psp . The relativistic correction 
is listed as 0{Za)^ in that reference. All further corrections 
to both the perturbative calculation and our calculation are 
contained in 'Remaining Corrections' which encompasses the 
muon AMM, proton MM, and all other lines of Table II in 
Ref. not already listed. We note the minor error in the last 
digit of the summary contribution of Ref. 0] for this state, 
likely arising from round-off error. 



Contribution 


Martynenko 


Present Work 


Dirac {V = Vc) 
Dirac {V = Vc + Vyp) 
Dirac {V = Vpc + Vvp) 




7.6141 
7.6159 
7.6204 


Dirac {V = Vpc + Vpvp) 
Leading contribution 
0{Za)'' contribution 


7.6018 
0.0011 


7.6204 


Subtotal 


7.6029 


7.6204 


Remaining Corrections 




0.3615 


Total: 


7.9644 


7.9819 



2P[=^ levels 



The 2P3/2 hyperfine structure, as derived 
in Ref. Q is given by 



Eq. (|4T|) corresponds to the leading term of Eq. (|36|) . to 
which the anomalous magnetic moments provide addi- 
tional corrections. Using the converged Dirac wave func- 
tions with exponential finite- Coulomb and finite vacuum 
polarization potentials (rather than Schrodinger wave 
functions) we calculate the expectation value of and 
find 



7.6204(5) meV. 



(42) 



The results of this calculation are summarised in Table lTvl 
where we note that the addition of the (point) vacuum 
polarization potential to the point-Coulomb potential in- 
creases the splitting by 0.0017(5) meV, and the introduc- 
tion of the finite- Coulomb potential increases this further 
by 0.0045(5) meV to arrive at the value above. The effect 
of finite- vacuum polarization is essentially zero here. 



C. 2P3/2 Hyperfine Structure 



2 

15 



30 



12mp(l -f- k) 



(43) 



where Ef is given in Eq. ((37)) (we note again the typo- 
graphical correction detailed in Section FVIII B|) . Alter- 
natively, inserting the relevant values for this state into 
Eq. (j39p and using the converged Dirac wave functions 
we find 



9 p 

AE^pl/' = 3.0415(5) meV 



(44) 



when the potential consists of the exponential finite- 
Coulomb and point vacuum polarization potentials. For 
this state, the addition of the (point) vacuum polariza- 
tion potential to the point-Coulomb potential increases 
the splitting by 0.0007(5) meV to the value listed in Ta- 
ble IYI and the introduction of the finite-Coulomb po- 
tential was found to have no effect within the limits of 
our calculation, so too was the introduction of the finite 
vacuum polarization potential. 



Following the same method as in the previous subsec- 

dF=: 
3/2 



tion, we can calculate the energy levels for the 2P^,7^^ and 
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TABLE V: Contributions to the 2P3/2 hyperfine splitting with 
comparison to values found in Martynenko 5] . Subscripts are 
defined in Table [T] Values are all in meV. Errors in the Dirac 
calculations are taken to be ±500 neV as per Section IIVI 
The details of the listed corrections are the same as those of 
Table IWl 



Contribution 


Martynenko 


Present Work 


Dirac {V = Vc) 
Dirac {V = Vc + VVp) 
Dirac {V = Vpc + Vyr) 




3.0408 
3.0415 
3.0415 


Dirac {V = Vpc + Vpvp) 
Leading contribution 
Relativistic correction 


3.0407 
0.0001 


3.0415 


Subtotal 


3.0408 


3.0415 


Remaining Corrections 




0.3518 


Total: 


3.3926 


3.3933 



IX. SUMMARY 

We summarise the findings of these non-perturbative 
Dirac calculations and compare to the previous literature 
values of perturbative calculations in Table I VII We add 
to this a combined expression for the cubic which when 
set equal to the experimental value of the measured tran- 
sition is solved to predict the proton rms charge-radius, 
as is done in Ref. 0- 

We do not include the result of hyperfine splitting cal- 
culations for the 2Pi/2 eigenstate as this is of no relevance 
to the measured transition. We also note the omission 
here of the energy shift attributed to a mixing between 
the 2P1/2 and 2P3/2 F = 1 states, as discussed in Ref. Q 
for comparison to Ref. where it is also absent. 



We find that the perturbative calculations are largely 
reproduced using our methods when considering the ap- 
propriate potentials for comparison. We further find that 
in several cases the use of the finite-vacuum polarization 
potential produces effects which are not accounted for in 
previous studies. The largest of these is the finite-size 
contribution to the 2Si/2 hyperfine splitting which has 
been neglected in the literature up to this point. 

Overall, the non-perturbative calculations do not elu- 
cidate any missing contributions of a magnitude large 
enough to resolve the proton radius problem outlined in 
Pohl et al. 0. 



X. CONCLUSIONS 

After careful consideration of the various contributions 
to the measured transition energy of Pohl et al. Q , cal- 



culated consistently using the Dirac equation with ap- 
propriate potentials, and following the addition of the 
required corrections to these calculations (taking further 
care to avoid overcounting issues), we find no single term 
which leads to a discrepancy with the perturbative re- 
sults of sufficient magnitude to account for the discrep- 
ancy reported in Ref. Q. These calculations nonetheless 
provide a useful insight into the reliability of the per- 
turbative calculations, and allow a simpler approach to 
future investigations. 

While it remains possible in principle that one or more 
of the higher-order corrections to the terms calculated in 
this work might be of sufficient magnitude to affect the 
analysis of Ref. Q, the precision with which the Dirac 
and perturbative calculations agree for the terms which 
we have calculated here strongly suggests that this will 
not be the case. 

In addition to the calculations presented here, we fur- 
ther note that our calculations of a hitherto overlooked 
contribution arising from off-mass-shell effects for the 
proton (which are negligible for electronic hydrogen) pro- 
vide a natural solution to the proton radius problem [l3| , 
and as such the combination of these two sets of cal- 
culations may be seen as a complete description of the 
measured transition in muonic hydrogen with no discrep- 
ancy in the rms charge radius of the proton. Because of 
the uncertain magnitude of the off-mass-shell effects, it 
is incorrect to complete the analysis of this transition to 
predict a proton rms charge radius — we await the results 
of current and future experiments which will be ascertain 
the strength of this contribution, after which a complete 
analysis will be possible. 

Nonetheless, we note that our calculations predict that 



the transition energy for the 2P^~^ 



to 2S 



F=l 
1/2 



transi- 



tion in muonic hydrogen is larger in magnitude than that 
which is predicted by the perturbative calculations, and 
that analysis of this data under the assumption that no 
further terms are required leads to the following values 
for the proton rms charge radius when fit to the experi- 
mental data; 



Pohl et al. 



Present Work 




0.84183(67) fm. 



0.83811(67) fm. 



The value listed as Present Work is taken as the solution 
to the cubic equation 

209.9505 - 5.2345(r2) + 0.0361 (r2)3/2 = 206.2949, (45) 

where the right-hand-side corresponds to the quoted 
value of the measured transition in Ref. Q ; the left- hand- 
side is taken from the relevant conclusion line of Table IVll 
and for which the errors in this calculation are domi nated 
by the experimental error. The extracted y^(r^) value 

listed as Pohl et al. is taken from Ref. (calculated 
in the same fashion) and differs from the central value 
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TABLE VI: Sum of perturbative and non-perturbative theoretical contributions to the measured experimental transition energy 
shown in Fig. [1] Subscripts are defined in Table [J Values are all in meV. The individual perturbative contributions (listed 
under 'Various') are taken fi:om Tables [iHVl In each case, the value given for the Dirac calculation is calculated using the 
combination of finite-Coulomb and finite vacuum polarization potentials {V — Vfc + Vfvp). The fractional factors for the 
hyperfine splittings are inserted for relevance to the measured transition, and are calculated via angular-momentum splitting 
rules. 



Contribution 


Various 


Present Work 


25'i/2-2Pi/2 Lamb shift (constant) 


206.0573 


206.0604 


2'S'i/2"2f'i/2 Lamb shift (finite-size) 


-5.2262(r2) -|- 0.0347(r|)3/2 


-5.2794(r2) -|- 0.0546(r2)3/2 


2P Fine Structure 


8.3521 


8.3521 


^ X 2Si^2 Hyperfine (constant) 


-5.7037 


-5.7345 


J X 2Si^2 Hyperfine (finite-size) 


0.0000 


0.0449 (r|) - 0.0185 (r-2)3/2 


1 X 2P3/2 Hyperfine 


1.2722 


1.2725 


(Various) Total (Perturbative) 
(Present Work) Total (Dirac) 


209.9779 - 5.2262 (: 
209.9505 - 5.2345 ( 


r^) -1- 0.0347 (r^)^/2 
r2) + 0.0361 (r^)3/2 
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quoted in 0; 0.84184(67), though the difference is weU 
within the quoted errors. 

We note the degree to which the cubic expression 
Eq. (|^5)) agrees with that of Ref. Q, despite the latter 
not involving a calculation of a finite-size contribution to 
the 25* hyperfine splitting. Some research in progress by 
the authors will elucidate some further overlooked contri- 
butions that will likely alter the agreement between these 
two expressions, and we look forward to future measure- 
ments with which we may compare our findings. 
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